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Abstract 

The GPS Master Control Station (MCS) currently makes significant use of the Allan Variance. 
This two-sample variance equation has proven excellent as a handy, understandable tool, both for 
time domain analysis of GPS Cesium frequency standards, and for fine tuning the MCS’s state 
estimation of these atomic clocks. 

The Allan Variance does not explicitly converge for the noise types of a < -3, and can be greatly 
affected by frequency drift. Because GPS Rubidium frequency standards exhibit non-trivial aging 
and aging noise characteristics, the basic Allan Variance analysis must be augmented in order to a) 
compensate for a dynamic frequency drift, and b) characterize two additional noise types, specifically 
a = -3 and a = -4. As the GPS program progresses, we will utilize a larger percentage of Rubidium 
frequency standards than ever before. Hence, GPS Rubidium clock characterization will require 
more attention than ever before. 

The three-sample variance, commonly referred to as a renormalized Hadamard Variance, is 
unaffected by linear frequency drift, converges for a > -5, and thus has utility for modeling noise in 
GPS Rubidium frequency standards. This paper demonstrates the potential of Hadamard Variance 
analysis in GPS operations, and presents an equation that relates the Hadamard Variance to the 
MCS’s Kalman Filter process noises (qs). 


INTRODUCTION 

The two-sample variance, or what we commonly refer to as the Allan Variance, has been an excellent 
device for time domain characterization of GPS Cesium frequency standards over the past few years. Over 
the past year, the GPS Master Control Station (MCS) has also applied the Allan Variance towards fine 
tuning the MCS’s state estimation of these Cesium clocks [3], 

In terms of Power-Law Spectral Density exponents, the Allan Variance does not explicitly converge for 
noise types of a < -3, and may be greatly affected by frequency drift [5], Because GPS Rubidium 
frequency standards exhibit significant aging and aging noise characteristics, the Allan Variance analysis 
must be augmented to dynamically compensate for this frequency drift, and to characterize two additional 
noise types, specifically a = -3 and a = -4. As the GPS program progresses, we will utilize a larger 
percentage of Rubidium frequency standards than ever before. In particular, the Block IIR satellite 
platform will house three atomic frequency standards, and two of these three will be Rubidium. Clearly, 
the characterization of GPS Rubidium clocks will soon require more attention than ever before. 
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In contrast, the three-sample variance, commonly referred to as a renormalized Hadamard Variance, is 
unaffected by linear frequency drift, converges for a > -5 [8], and hence has a potential utility for modeling 
the various noise types resident in GPS Rubidium frequency standards. This paper demonstrates this 
potential for Hadamard Variance analysis in GPS analysis operations, and presents the relationship 
between the Hadamard Variance and the MCS’s Kalman Filter process noises (q's). 


THE HADAMARD VARIANCE EQUATION 

A mainstay of atomic clock characterization, the two-sample (Allan) Variance essentially examines the 
second difference of phase, equivalent to the first difference of the time-averaged frequencies over two 
successive adjacent time intervals (r) [5]: 
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Similar in principle to the structure to the Allan Variance, the three-sample variance examines the third 
difference in phase, equivalent to the second difference of the time-averaged frequencies over three 
successive adjacent time intervals (r). The timing community has commonly referred to this three-sample 
variance as the Hadamard Variance. Though the term Hadamard Variance has been used more generally 
in various applications of multi-sample time domain analysis, for the purposes of this paper, we shall 
define the Hadamard Variance as follows: 


\ M ~ 2 — _ _ 2 _ 

= ^^](y i+2 -2y i+ , + y ( ) , y, = the time - averaged frequency over r ( . (2) 

The Hadamard Deviation (the square root of the Hadamard Variance) identifies two noise types that the 
Allan Deviation does not explicitly identify [8], For this paper, we shall name the following noise types: 
for a = -3, “Flicker Walk FM”; for a — - 4, “Random Run FM” [5], Figure 1 visually describes the noise 
types identified by the Hadamard Deviation [6], 

In terms of phase, equation (2) converts to [8]: 

| N-3 2 

„g%(r) = 2 - -- V(x )>3 -3x m +3*,+, -x,) , x, = the phase measurement at t r (3) 

or (TV — j) 


or, equivalently: 


H <j2 y( T ) = E [( X ,+ 3 -*, +2 )-(x i+2 -* <+I )-(x f+2 -*„) + (*„ -x,)] 2 (4) 

where E[.] is the expectation operator. Each phase measurement x, = x(t,) in equation (4) is separated from 
each neighboring successive phase measurement by a time interval value of r. Meaning, 


292 



x(t,+i) = x(t, + l), 

x(t, +2 ) = x(t, +! + x) = x(t, + 2i), and 
x(/, + i) = x(t i+2 + r)= x(t i+! + 2r) = x(f, + 3r) 


(5) 


( 6 ) 

(7) 


MCS KALMAN FILTER TIME UPDATE PREDICTIONS 

The propagation (time update), of Rubidium clock states in the MCS Kalman Filter, is modeled using the 
following polynomial expansion [7]: 
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where t is the prediction span, and x(t), y(t), and z(t) are the phase, frequency, and frequency drift values, 
respectively, of the clock in question. Note that y(t) is the time derivative of x(t), and z(t) is the time 
derivative ofy(0 A(x), My), and A(z) are assumed to be random error increments, independent of x(t), y(t), 
and z(f), having a prediction covariance P represented by a function of the Kalman Filter process noises 
(^s) [1,7]: 
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An expansion of equation (8) produces the following equations: 
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Using equations (10), (1 1), and (12), and examining the differences between each successive x,: 
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(13) 


(x i+ , - xl) = ry,+ (1 IT)?Zi + Ax t +, 

(x i+2 - x 1+/ ) = ry t+ i + (1 12)?z i+ , + Ax H2 (14) 

(x i+3 - X I+2 ) = T y i+2 + (1 H)^Zi + 2 + A x i+3 (15) 

AN EXPECTATION OPERATOR EXPANSION OF ho^t) 

Inserting equations (13), (14), and (15) into the following expression: 

[ThirdDiJ] = [(x, +2 - x, +2 ) - (x /+2 - x, +; ) - (x l+2 - x 1+ ,) + (x, +; - x,)] (16) 

obtains: 

[ThirdDiJ] = Ax I+J - 2Ax, +2 + Ax <+/ + r{(y, -y i+l ) - (y l+ , -y, +2 )} + (1/2)i 2 {(z, - z i+I ) - (z, +/ - z, +2 )} (17) 

Examining the differences between each successive^/ and z„ from equations (10), (1 1), and (12): 

(yi+i ~yd = + Ay i+I (18) 

O', +2 ->'i + ^) = tZ/ + ; + A y i+ 2 (19) 

(Z/+; * Z|) = AZ/+/ (20) 

(z I+2 -z, +/ ) = Az 1+2 (21) 

Equation (17) translates into: 

[ ThirdDiJ] = Ax 1+2 - 2Ax l+2 + Ax, +i + t{{tz, +i + A y i+2 ) - (rz, + A y i+ ,)} + (l/2)t 2 {(Az/ +2 ) - (Az I+ ,)} (22) 

With some more algebraic manipulation: 

[■ ThirdDif] = Ax, +i - 2Ax (+2 + Ax, +; + r{(Ay, +2 - Ay i+/ ) + (tz, +/ - rz,)} + (l/2)r 2 {(Az i+2 ) - (Az 1+/ )} (23) 

[ThirdDiJ] = Ax, +J - 2Ax (+2 + Ax j+/ + r{(A>>, +2 - Ay I+/ ) + (rAz i+/ )} + (1/2 )i 2 {(Az <+ 2 ) - (Az (+/ )} (24) 

[ThirdDiJ] = Ax I+J - 2Ax, +2 + Ax, +; + r{(Ay/ +2 - Ay 1+ ,)} + (l/2)r 2 {(Az, +2 ) + (Az, + ,)} (25) 

[ThirdDiJ] = {Ax 1+2 } + {-2Ax/ +2 + rA y H2 + (l/2)t*(A z i+2 )} + {Ax I+ , - rAy i+l + (l/2)i*(Az,+/)} (26) 

Since we’ve now broken down this expansion into three independent polynomial terms: 

E[ThirdDiJ] 2 = E[{Ax, +J } + {-2Ax, +2 + rAy< +2 + (l/2)«*(Az <+2 )} + {Ax, +/ - rA y i+l + (l/2)i ? (Az< +/ )}] 2 (27) 

the independence of each term {} allows us to separate equation (27) into three individual expectation 
operators [4]: 
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E[ThirdDij f = E[Ax, +J ] 2 + E[-2Ax 1+? + rA y I+2 + (l/2)r(Az I+z )] 2 + E[Ax i+ , - rAy I+ , + (l/2)^(Az, +; )] 2 (28) 
Expressing each term ofE[ThirdDi/\ 2 as a function of the Kalman Filter prediction covariance matrix [1]: 
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By adding each term: 

E[ThirdDif] 2 = 6^i(r) + ^(z 3 ) + (1 1/20)^ 3 (t^) (32) 

Hence, 

H CF 2 y (i) = -^-rE[ThirdDiff = — ( r ) + (r 3 ) + (1 1 / 20)^ 3 (r 5 )] (for a = 0, -2, -4) [6] (33) 

6 t 6 r 

„o* y (T) = q x x A + (l/6)q 2 T+ (1 1/120)^3^ (for a = 0, -2, -4) (34) 


RELATING WHITE PM TO THE HADAMARD VARIANCE 

Equation (34) does not, however, account for white PM noise (a = 2) [6], sometimes also referred to as 
representation error [1]. The Hadamard Variance can be expressed, in terms of phase measurements, as 
follows: 


„cr%(r) = — -y E[x 1+3 -3x, + 2 + 3x l+1 -xj 2 (35) 

or 

When white PM is the only significant noise component, the individual phase values are uncorrelated with 
time, and may be separated [4]: 
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(36) 


H <7 2 y(r) = -L{e[x, +3 ] 2 + E[-3x, +2 ] 2 + E[3x I+1 ] 2 + E[-x,] 2 } (for a = 2) 

When assuming that white PM is the primary noise source, the representation error, which we’ll denote as 
q 0 = E[x,] 2 , is independent of t„ and thus is a phase variance that is constant across time. Therefore, 

h ° 2 r(T ) = {(1 + 9 + 9 + l)E[x, ] 2 } = j^q 0 (for a = 2) (37) 


THE HADAMARD-Q EQUATION 

In the presence of both a) white PM, and b) the three noise types modeled by P, and assuming independence 
between the white PM and the other noise types, equations (34) and (37) can be combined into one that 
models four noise types, namely a = 2, 0, -2, and -4: 

«0y(r) = (10/3)? 0 r 2 + q x f x + (l/6)? 2 r+ (1 1/120)^ (38) 

Note how this equation compares to the analogous equation relating the Allan Variance to the qs [2,3]: 

cfyir) = 3q 0 f 2 + q x r A + (l/3)^ 2 r+ (1/20)^^ (39) 

For white FM, the Allan and Hadamard Variances are mathematically equivalent. For white PM, the two 
Variances are roughly the same, and, for random walk FM, the Variances differ by a factor of two. 
Though an analyst may use either the Hadamard Variance or the Allan Variance for deriving MCS qs, each 
has its own set of advantages and disadvantages. 

The primary advantage of the Hadamard Variance is the automatic removal of linear frequency drift [8], 
Whereas the equation relating the Allan Variance to MCS qs assumes that the analyst must apply a 
continuously dynamic correction for frequency drift, the Hadamard-Q equation doesn’t require this 
assumption. The tradeoff, however, is that the Hadamard Variance incurs an extra computational burden, 
simply because it examines the third (vice the second) difference of phase. For analyzing GPS Cesium 
frequency standards, the increased computational load of the Hadamard Variance proves fruitless, only 
because the Allan Variance gets the job done more efficiently [3], 

Many timing experts, over the years, have extensively used techniques for applying a continuously dynamic 
correction for frequency drift, prior to using the Allan Variance for deriving q 3 values with high confidence. 
This paper does not address the issue of confidence in the qi value produced by the Hadamard-Q equation. 
This paper does, however, present an easily understood relationship between a relatively lesser known 
equation (the Hadamard Variance), and a set of system parameters used by the MCS (the Kalman Filter 
< 7 s). On initial appearance, given the computational capability, one can see the great potential utility of an 
algorithm that applies a relatively simple equation onto a large measurement data base, in order to derive 
Kalman Filter qs, without the need to apply preparatory frequency drift corrections. In the future, the 
author hopes to further investigate a) the real-world utility of this relationship, b) the issue of estimate 
confidence, and c) the net gain from the increased utility balanced against the increased computational 
burden. 
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CONCLUSION 


The implication of the Hadamard Variance in GPS operations is as follows: Analysts at the MCS could 
simply gather a large data base of clock phase measurements, apply all known step corrections, perform a 
number of iterations of the Hadamard Variance equation, and plot the results to visually describe the noise 
characteristics of GPS Rubidium clocks (including a = -4). Consequently, the Hadamard Variance could 
offer GPS operators an alternate tool for characterizing GPS atomic frequency standard noise. 

Perhaps more significantly, the implication of equation (38) is that GPS analysts now have an easily 
understood technique to relate raw clock phase measurements towards deriving important Kalman Filter 
clock estimation parameters ( qs ), that are unique to the performance of each individual clock, and that will 
include qi automatically, without any need for the preliminary removal of frequency drift. 

The MCS hopes to make continued use of the Allan Variance, for both the characterization of Cesium 
clocks, and the derivation of their associated process noise values. This application of the Hadamard 
Variance widens the array of available tools for the characterization of Rubidium clocks, and the derivation 
of their associated process noise values. The GPS Block IIR satellite program will use a large percentage 
of Rubidium clocks. Certainly, the ever-important issue of refining Rubidium clock estimation may see its 
most important days in the years ahead. 
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Noise Types Identified by the Hadamard Deviation 



Noise Type Power-Law Spectral Density Exponent hOvCt) Log-Log Slope 


White PM 

a = 2 

-i 

Flicker PM 

a = 1 

-i 

White FM 

a = 0 

- 1/2 

Flicker FM 

a = -1 

0 

Random Walk FM 

a = -2 

1/2 

Flicker Walk FM 

a = -3 

1 

Random Run FM 

a = -4 

3/2 


Figure 1 




Questions and Answers 

SERGEY V. ERMOLIN (HEWLETT-PACKARD): The Hadamard Variance does remove 
linear drift, that’s true; and that saves you some time on preprocessing. But, it doesn’t remove 
any drift beyond linear. As you showed on one of your first slides, that rubidium standards on 
board one of the space vehicles show not only linear drift, but possibly quadratic drift to some 
other power. 

So still, if you wanted to go ahead with longer operational time, you would still have to do 
some preprocessing, even if you used the Hadamard Variance. 

CAPTAIN STEVEN HUTSELL (USAF): The analogy is the Allan Variance does not 
care about a constant frequency offset. For instance, our atomic clocks can have a 1 x 10 -11 
frequency offset; but if it’s stable enough, the Allan Variance will be low, regardless of that 
frequency offset. The analogy is in the Hadamard Variance, if the satellite clock, or whatever, 
has an already existing common offset of, say, 3 x 10 -18 seconds per second squared of frequency 
drift, that will not adversely affect the Hadamard Variance calculations. 

However, in the Allan Variance, if there is a frequency drift, it will affect it. But with the 
Hadamard, it won’t. In the same sense, a random walk in random walk FM will affect the 
Allan Variance by causing a positive one slope. But the frequency offset itself will not affect it. 

DAVID ALLAN (ALLAN’S TIME): Dr. Barnes did some work three decades ago on the 
confidence question, actually related to it, in the third difference estimate. The bottom line 
was that the confidence is worse by a significant amount, especially when you have finite data 
lengths that come into impact you quite adversely. 

The other point is that it turns out a logarithmic drift estimator, both for quartz and rubidium 
is much better than linear. I think the graph that you showed outside of the turn-on transient 
probably would fit a logarithmic curve quite well. 

So, one might be better doing logarithmic modeling if we deal with a lot of rubidiums in the 
future. Yes, one would expect that the logarithmic function could be fit to this quite well, 
outside of the first point. 

The last one I would like to suggest that you think about — and I want to say that I think 
you’ve done a beautiful piece of work, but you can actually remove the effect of drift, kind of 
in real time, from the second difference operator, because you know the exact equation for 
the effect of drift on it until you can subtract that from and get an estimator variance without 
the drift effectively in real time. So, it doesn’t need to impact the value of the variance if you 
don’t want it to; and it gives you a tighter confidence of the estimate. 

CAPTAIN STEVEN HUTSELL (USAF): Yes, and the intent of this is not really to present 
the best way to estimate frequency drift dynamically. Really, what I wanted to do was examine 
the way the MCS is currently set up. And right now, it does not have the capability to do 
what you just described. The MCS is only set up as a dynamic Kalman filter three-state vector 
that needs process noise values. 

I completely agree that are far more sophisticated techniques to look at frequency drift than 
what’s set up in the Kalman filter. Sorry, your first comment? 
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DAVID ALLAN (ALLAN’S TIME): [Inaudible]. 

CAPTAIN STEVEN HUTSELL (USAF): Yes, we see that too. We see it start to converge 


DAVID ALLAN (ALLAN’S TIME): [Inaudible]. 

CAPTAIN STEVEN HUTSELL (USAF): We agree. Over time, the frequency drift goes 
from a negative value, around 3 or 4 x 1CT 18 , and gradually starts logarithmically to approach 

zero. 

At the beginning, however, sometimes we see it start hugely negative, like —1 x 10 1 . Sometimes 

+1 x 1(T 1G . We’re talking about over the first 48 hours that we turn it on. We would need to 
address how we try to model that. It’s also probably appropriate to ask what’s causing that, 
what physical phenomenon is causing it to be positive at the beginning for some clocks, and 
negative for the others. But, I do agree. 
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